#ifndef _CFD_UTILS_H
#define _CFD_UTILS_H

/*
 * The following file contains all the extra macros/data structures/functions
 */

/********************************************************
 * 			Macros
 *******************************************************/

#define 	TOL 		10e-2	/**< Tolerance used for comparisons */


#define		NO_SLIP		1
#define		FREE_SLIP	2
#define		OUTFLOW		3

#define		C_B		0x00000100	/**< The obstacle cell */
#define 	C_F		0x00000200	/**< The fluid cell */
#define		C_Y		0x00000300	/**< The Empty cell */

#define		MASK		0x00000300	/**< Mask used in check the cell type */

#define		I_MASK		0x000000FF	/**< Mask used to determine interior cells */

#define		CLEAR		0x00000F00	/**< Mask used to clear the neighbors in mark_cells */

/*
 * The boundary markers
 *
 * In order to be able to support three different cases (obstacle, fluid, and empty),
 * each of the fields has been extended to be 2-bits in length.
 *
 * The following scheme is used
 * 00 = Not Used
 * 01 = Obstacle
 * 10 = Fluid
 * 11 = Empty
 *
 * So for each of the fields, we just plug-in the value corresponding to its condition.
 *  ----------------------------------------------------------------
 *  | Left P | Right P| Center |  East  |  West  | South  | North  |
 *  ----------------------------------------------------------------
 *  | 1-bit  | 1-bit  | 2-bits | 2-bits | 2-bits | 2-bits | 2-bits |
 *  ----------------------------------------------------------------
 *
 * All the macros' values are given in hexadecimal to make it simpler for
 * debugging mis-placed bits.
 */
#define 	B_N		0x00000102	/**< North edge cell */
#define 	B_S		0x00000108	/**< South edge cell */
#define		B_W		0x00000120	/**< West  edge cell */
#define		B_E		0x00000180	/**< East  edge cell */

#define		B_NE		0x00000182	/**< North-East corner cell */
#define		B_NW		0x00000122	/**< North-West corner cell */
#define		B_SE		0x00000188	/**< South-East corner cell */
#define		B_SW		0x00000128	/**< South-West corner cell */

/* The masks used in setting the neighbors */
#define 	F_N		0x00000002	/**< Fluid in the North */
#define 	F_S		0x00000008	/**< Fluid in the South */
#define		F_W		0x00000020	/**< Fluid in the West */
#define		F_E		0x00000080	/**< Fluid in the East */

#define		Y_N		0x00000003	/**< Air in the North */
#define		Y_S		0x0000000C	/**< Air in the North */
#define		Y_W		0x00000030	/**< Air in the North */
#define		Y_E		0x000000C0	/**< Air in the North */


/*
 * The empty neighbours markers
 */
#define		C_N		0x00000203
#define		C_S		0x0000020C
#define		C_W		0x00000230
#define		C_E		0x000002C0

#define		C_SE		0x000002CC
#define		C_SW		0x0000023C
#define		C_NS		0x0000020F
#define		C_WE		0x000002F0
#define		C_NE		0x000002C3
#define		C_NW		0x00000233

#define		C_NSE		0x000002CF
#define		C_NSW		0x0000023F
#define		C_NWE		0x000002F3
#define		C_SWE		0x000002FC

#define		C_NSWE		0x000002FF


/* The left and right pressure flags */
#define		P_L		0x00000800	/**< The left boundary pressure flag */
#define		P_R		0x00000400	/**< The right boundary pressure flag */

/*
 * Problem definitions
 */
#define		KARMAN 			"karman"	/**< Von Karman Vortex Street Problem */
#define		CAVITY			"cavity"	/**< Driven Cavity Problem */
#define		PLANE			"plane"		/**< Plane Shear flow Problem */
#define		STEP			"step"		/**< Step Problem */
#define		NOZZLE			"nozzle"	/**< Nozzle Problem */
#define		REVERSE_FLOW		"reverse_flow"	/**< Reverse Flow Problem */
#define		FALLING_DROPS		"Falling_Drops"	/**< Falling Drops Problem */
#define		BREAKING_DAM		"Breaking_Dam"	/**< Breaking Dam Problem */
#define		PROBLEM_MAX_LENGTH	16		/**< The maximum length for the problem name */

/*
 * Flow visualization parameters
 */
#define		STREAK_LINES_PART_MAX	1000		/**< The maximum number of particles visualized using the streaklines */


/*************************************************************
 * 			Data Structures
 ************************************************************/

/**
 * Holds all the data related to the geometry info.
 */
struct geometry_data
{
	double xlength; 	/**< Domain size in x-direction */
	double ylength; 	/**< Domain size in y-direction */
	int imax;  		/**< Number of interior cells in x-direction */
	int jmax; 		/**< Number of interior cells in y-direction */
	double dx; 		/**< length dx of one cell in x-direction */
	double dy; 		/**< length dy of one cell in y-direction */
};

/**
 * Holds all the info related to time-stepping and simulation duration
 */
struct time_data
{
	double t; 		/**< The current time value */
	double t_end; 		/**< The final time */
	double dt;  		/**< The time step size */
	double tau; 		/**< The safety factor for the time step control */
	double dt_value; 	/**< The time interval for writing the visualization data into the output file */
};

/**
 * Holds the data related to solving the pressure equation using the SOR solver
 */
struct pressure_iteration
{
	int itermax;		/**< The maximum number of pressure iterations in one time step */
	int it;  		/**< SOR iteration counter */
	double res; 		/**< The residual norm of the pressure equation */
	double eps; 		/**< The accuracy criterion (tolerance) for pressure equation (res < eps) */
	double omg;		/**< The relaxation factor for the SOR solver */
	double alpha;		/**< The upwind differecing factor */
};

/**
 * Holds the problem dependent quantities
 */
struct problem_parameters
{
	double Re;		/**< Reynolds number */
	double GX; 		/**< X-component of the external forces  */
	double GY;		/**< Y-component of the external forces  */
	double UI;	 	/**< Initial data for U */
	double VI;	 	/**< Initial data for V */
	double PI;	 	/**< Initial data for P */
};

/**
 * Holds the boundary options switch
 */
struct boundary
{
	int wl;					/**< The left boundary switch */
	int wr; 				/**< The right boundary switch */
	int wt;					/**< The top boundary switch */
	int wb;	 				/**< The bottom boundary switch */
	double dp;				/**< The pressure difference */
	char problem[PROBLEM_MAX_LENGTH];	/**< The problem name */
	double driving_velocity;		/**< The value of the velocity */
	int ppc;				/**< The number of particles per cell in each dimension */
	int N;					/**< Number of particles lists */
	int npart;				/**< The total number of fluid particles */
};

/**
 * Holds the parameters used for the flow visualization (pathlines and streaklines)
 */
struct flow_vis
{
	int N;			/**< The number of particles equidistantly distributed over the line specified by (x1,y1) and (x2,y2) */
	double x1;		/**< The x-coordinate of the 1st point in the line used for particles/streaklines tracing */
	double y1;		/**< The y-coordinate of the 1st point in the line used for particles/streaklines tracing */
	double x2;		/**< The x-coordinate of the 2nd point in the line used for particles/streaklines tracing */
	double y2;		/**< The y-coordinate of the 2nd point in the line used for particles/streaklines tracings */
	double dt_path;		/**< The time step used to output path-lines */
	double dt_streak;	/**< The time step used to output streak-lines */
	double dt_insert;	/**< The time step used to insert new particles in the case of streak-lines */
};

double **U; 			/**< The velocity in the x-direction */
double **V; 			/**< The velocity in the y-direction */
double **P; 			/**< The pressure */
double **RS; 			/**< The right-hand side of the pressure equation */
double **F; 			/**< F */
double **G; 			/**< G */

int **Flag; 			/**< Flag */

int nCount;			/**< Used for particles tracing */




/**************************************************************
			Functions
**************************************************************/

/**
 * Checks whether a given cell is a fluid cell or not
 * @param	flag	The value of Flag[i][j]
 * @return	1 if the cell represents a fluid cell, 0 otherwise
 */
inline int is_fluid(int flag);


/**
 * Checks whether a given cell is an obstacle cell or not
 * @param	flag	The value of Flag[i][j]
 * @return	1 if the cell represents an obstacle cell, 0 otherwise
 */
inline int is_obstacle(int flag);

/**
 * Checks whether a given cell is an empty cell or not
 * @param	flag	The value of Flag[i][j]
 * @return	1 if the cell represents an empty cell, 0 otherwise
 */
inline int is_empty(int flag);

/**
 * Checks whether a given cell is an interior cell or not
 * @param	flag	The value of Flag[i][j]
 * @return	1 if the cell represents an empty cell, 0 otherwise
 */
inline int is_interior(int flag);

#endif
